close all;clear all;clc;
% Parameters
A = 1.1;
g = 2;
lambda = 0.2;
eta = 0.6;
gprime = 1;
k = 0.5;

% Functions
tauf = @(g, gprime, lambda, eta, A, k)...
    1/(g-gprime)*(log(A)-log((exp((g-gprime)*k*eta)-1)/((g-gprime)*eta)...
    + (1-k)*((g/lambda)*exp((g-lambda-gprime)*k*eta)...
    - (g-lambda)/lambda*exp((g-gprime)*k*eta))));

alphaf = @(g, gprime, eta, A, k, tau)...
    (A*exp(gprime*(tau+k*eta))...
    - exp(g*tau)*(exp(g*k*eta)-exp(gprime*k*eta))/((g-gprime)*eta))/(1-k);

hf = @(g, lambda, eta, k) lambda*exp(lambda*k*eta)/(exp(lambda*k*eta)-1);

kcf = @(g, gprime, lambda, eta) fsolve(@(k) (1-k)*g*eta + g/lambda...
    - (g-lambda)/lambda*exp(lambda*k*eta) - exp(-(g-lambda-gprime)*k*eta), 0);

% Calculations
tau = 0:0.01:0.4;

for i = 1:length(tau)
    alpha(i) = alphaf(g, gprime, eta, A, k, tau(i));
    h(i) = hf(g, lambda, eta, k);
    decomp1_h(i) = 1/g*log(h(i)/(h(i)-g));
    decomp3_alpha(i) = 1/g*log(alpha(i)*exp(-g*k*eta));
    LHS(i) = tau(i);
    RHS(i) = decomp1_h(i) + decomp3_alpha(i);
end

% points for figures
for j = 1:length(tau)
    if (j==1)|(j==length(k))|(mod(j,4)==0)
        tau_figure(j) = tau(j);
        decomp3_alpha_figure(j) = decomp3_alpha(j);
        decomp1_h_figure(j) = decomp1_h(j);
    else
        tau_figure(j) = NaN;
        decomp3_alpha_figure(j) = NaN;
        decomp1_h_figure(j) = NaN;
    end
end

knc = min(1/(lambda*eta)*log(g/(g-lambda)),...
    1/((g-gprime)*eta)*log(A*(g-gprime)*eta+1));

taustar = tauf(g, gprime, lambda, eta, A, k);


% Figures
h = figure;
set(h,'Position',[0 0 800 800]);
set(h, 'defaultLegendInterpreter','latex',...
    'defaultAxesTickLabelInterpreter','latex');
a = axes;

p1 = plot(tau, decomp3_alpha_figure, '-.k*', 'linewidth', 2, 'MarkerSize', 8);
hold on;
p2 = plot(tau, decomp1_h_figure, '-.ok', 'linewidth', 2, 'MarkerSize', 8);
p3 = plot(tau, RHS, 'k', 'linewidth', 2);
p4 = plot(tau, LHS, 'k--', 'linewidth', 2);
plot(tau, decomp3_alpha,'-.k','linewidth',1);
plot(tau, decomp1_h,'-.k','linewidth',1);
xl = xlim();
yl = ylim();
plot([taustar taustar], [yl(1) taustar], 'k--', 'linewidth', 1);
plot([xl(1) taustar], [taustar taustar], 'k--', 'linewidth', 1);
plot([xl(1) xl(2)], [0 0], 'k--', 'linewidth', 1);

set(a, 'YTick',[-0.4 0 taustar 0.4],...
    'YTicklabel',{'$-0.4$' '$0$' '$\tau^{*}$' '$0.4$'},...
    'XTick',[0 taustar 0.4],...
    'XTicklabel',{'$0$' '$\tau^{*}$' '$0.4$'}, 'Fontsize', 16);

text(0.07, 0.03, '$45^{\circ}$', 'interpreter','latex', 'Fontsize', 20);

xlabel('Waiting Time ($\tau$)', 'interpreter','latex');
legend1 = legend([p1 p2 p3],...
    {'$\frac{1}{g}\log\left[\alpha(\tau,k) e^{-g k \eta}\right]$',...
    '$\frac{1}{g}\log\left(\frac{h(k)}{h(k)-g}\right)$',...
    '$\frac{1}{g}\log\left[\alpha(\tau,k) e^{-g k \eta}\right]+\frac{1}{g}\log\left(\frac{h(k)}{h(k)-g}\right)$'},...
    'location', 'southwest');

%saveas(h, 'tau_decomposition_tau.eps');



